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Abstract  We  provide  an  analytical  study  on  the  stability  of  equilibria  of  rigid  rodlike 
nematic  liquid  crystalline  polymers  (LCPs)  governed  by  the  Smoluchowski  equation  with 
the  Maier-Saupe  intermolecular  potential.  We  simplify  the  expression  of  the  free  energy 
of  an  orientational  distribution  function  of  rodlike  LCP  molecules  by  properly  selecting 
a  coordinate  system  and  then  investigate  its  stability  with  respect  to  perturbations  of 
orientational  probability  density.  By  computing  the  Hessian  matrix  explicitly,  we  are 
able  to  prove  the  hysteresis  phenomenon  of  nematic  LCPs:  when  the  normalized  polymer 
concentration  b  is  below  a  critical  value  b*  (6.7314863965),  the  only  equilibrium  state  is 
isotropic  and  it  is  stable;  when  b*  <  b  <  15/2,  two  anisotropic  (prolate)  equilibrium  states 
occur  together  with  a  stable  isotropic  equilibrium  state.  Here  the  more  aligned  prolate 
state  is  stable  whereas  the  less  aligned  prolate  state  is  unstable.  When  b  >  15/2,  there  are 
three  equilibrium  states:  a  stable  prolate  state,  an  unstable  isotropic  state  and  an  unstable 
oblate  state. 
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1  Introduction 

Since  its  discovery  by  Austrian  botanical  physiologist  Friedrich  Reinitzer  in  1888,  liquid 
crystals  have  spurted  intensive  experimental,  theoretical  and  numerical  studies  [1,  5,  6,  18,  31, 
39].  One  notable  example  is  that  in  1991  Pierre-Gilles  de  Gennes  was  awarded  the  Nobel  Prize 
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in  Physics  “for  discovering  that  methods  developed  for  studying  order  phenomena  in  simple 
systems  can  be  generalized  to  more  complex  forms  of  matter,  in  particular  to  liquid  crystals 
and  polymers”. 

Liquid  crystals  refer  to  a  state  of  matter  that  has  properties  between  those  of  a  solid  crystal 
and  those  of  an  isotropic  liquid.  For  example,  a  liquid  crystal  may  flow  like  an  ordinary  liquid, 
but  its  molecules  may  be  oriented  in  some  directions.  Examples  of  liquid  crystals  can  be  found 
both  in  nature  (e.g.  solutions  of  soap  and  detergents)  and  in  technological  applications  (e.g. 
electronic  displays) .  There  are  many  different  types  of  liquid  crystal  phases,  among  which  one 
of  the  most  common  phases  is  the  nematic.  In  a  nematic  phase,  the  rod-shaped  molecules  have 
no  positional  order  but  they  display  long-range  directional  order  with  their  long  axes  roughly 
aligned  parallel  to  a  common  axis  (which  is  called  a  “director”).  Most  nematics  are  uniaxial, 
meaning  that  they  have  one  axis  that  is  longer  and  preferred  while  the  other  two  axes  are 
equivalent.  However,  some  nematics  are  biaxial  which  implies  that  the  molecules  also  orient 
along  a  secondary  axis  in  addition  to  orienting  their  long  axis.  Nematics  can  be  easily  aligned 
by  applying  an  external  magnetic  or  electric  field.  The  optical  properties  of  aligned  nematics 
make  them  extremely  useful  in  liquid  crystal  displays. 

Liquid  crystallinity  in  polymers  may  occur  either  by  heating  a  polymer  above  its  melting 
point  (thermotropic  LCPs)  or  by  dissolving  a  polymer  in  a  solvent  (lyotropic  LCPs).  Liquid 
crystalline  polymers  combine  the  properties  of  liquid  crystals  with  the  properites  of  polymers. 
As  a  result,  it  is  easy  to  guide  LCPs  by  applying  small  external  forces  such  as  flows  and  it  is 
easy  to  control  their  chemical  and  material  parameters.  LCPs  exhibit  certain  anisotropy  of  the 
electrical,  mechanical  and  magnetic  properties  with  high  flexibility  or  elasticity.  They  are  also 
exceptionally  unreactive  and  inert,  and  highly  resistant  to  fire.  Due  to  their  various  properties, 
LCPs  have  found  more  and  more  applications,  including  electrical  and  mechanical  parts,  food 
containers,  and  any  other  applications  requiring  chemical  inertness  and  super  strength.  A  main 
example  of  LCPs  in  solid  form  is  high  strength  fibers  which  have  been  made  into  bullet-proof 
vests  and  airbags. 

The  theoretical  studies  of  liquid  crystals  traced  back  more  than  60  years  ago.  In  1949 
Onsager  [29]  developed  a  statistical  theory  showing  that  simple  exclude-volume  repulsions  be¬ 
tween  long  rods  are  sufficient  to  create  a  liquid  crystal.  Onsager’s  theory  is  valid  for  dilute 
liquid  crystals.  Flory  [9,  10]  presented  a  lattice  model  and  his  theory  yields  good  approxima¬ 
tion  for  dense  and  highly  ordered  liquid  crystals.  Both  the  Onsager  and  Flory  theories  predict 
an  isotropic-nematic  phase  transition  when  the  concentration  is  high  enough.  In  three  pub¬ 
lications  Maier  and  Saupe  [26-28]  proposed  a  mean-field  theory  to  explain  the  occurrence  of 
a  liquid  crystal  phase  in  low  molecular  weight  materials  based  on  the  assumption  that  only 
induced  dipolar  forces  are  relevant  for  orientational  ordering.  An  impressive  signature  for  the 
popularity  of  the  Maier-Saupe  theory  is  its  high  number  of  citations  (almost  4000  up  to  now) 
over  the  years  [30].  In  1965  Landau  and  Ginzburg  [22]  developed  a  phenomenological  theory  for 
liquid  crystals.  Even  though  it  has  success  in  most  cases,  Landau-Ginzburg  theory  is  not  exact 
and  it  suffers  failure  in  some  cases  [18].  One  of  the  earliest  constitutive  theories  of  liquid  crystals 
is  the  Leslie-Ericksen  vector  theory  [7,  23]  in  1970s.  The  Leslie-Ericksen  theory  is  attractive  for 
studying  low  molar-mass  liquid  crystals. 

After  Maier  and  Saupe  and  Flory’s  work,  many  new  theories  have  been  developed  for  liquid 
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crystalline  polymers.  One  of  the  most  famous  kinetic  models  for  polymer  liquid  crystals  was 
developed  by  Doi  and  Edwards  [6]  in  1980s  and  Hess  [19].  The  crucial  idea  of  the  Doi  theory  is  to 
model  the  polymeric  molecules  as  a  suspension  of  rigid  rodlike  nematogens  and  then  describe  the 
ensemble  with  an  orientational  probability  distribution  function.  The  orientational  distribution 
function  is  governed  by  a  nonlinear  Smoluchowski  (or  Fokker-Planck)  equation  which  takes  into 
account  the  hydrodynamic,  Brownian  and  intermolecular  forces. 

In  recent  years,  the  Doi  theory  has  attracted  enormous  mathematical  attentions  (for  ex¬ 
ample,  [2-4,  8,  11-17,  21,  24,  25,  32-38]).  One  mathematical  advance  in  the  understanding 
of  the  equilibrium  states  of  the  Doi  model  is  the  rigorous  proof  on  that  the  equilibria  of  the 
Smoluchowski  equation  with  the  Maier-Saupe  potential  are  uniaxial  [8,  24,  35].  However,  it 
seems  more  challenging  to  obtain  a  solid  proof  on  the  stability  of  the  equilibria.  The  purpose 
of  this  paper  is  to  circumvent  this  difficulty  and  seek  a  mathematical  proof  on  the  stability  of 
the  equilibria.  The  main  result  proven  here  is:  all  stable  equilibria  are  either  isotropic  or  highly 
aligned  prolate  uniaxial. 


2  Background  and  Formulation 


For  reader’s  convenience,  we  start  with  the  briefest  glimpse  of  the  Doi  kinetic  theory. 

Let  us  denote  the  orientational  direction  of  each  rodlike  molecule  by  a  unit  vector  m, 
the  probability  density  function  by  p(m).  The  evolution  of  the  probability  density  function  is 
described  by  the  Smoluchowski  equation  [6] : 


dp  d  (  1  dU  dp  \ 
dt  dm  \ksT  dmP  +  dm)  ' 


(1) 


where  d/dm  is  the  orientational  gradient  operator  [1],  D  the  rotational  diffusivity,  kB  the 
Boltzmann  constant,  and  T  the  absolute  temperature.  For  simplicity,  we  set  kBT  =  1  which  is 
equivalent  to  normalizing  U  by  kBT. 

For  nematic  polymers,  the  total  potential  consists  of  only  the  Maier-Saupe  interaction 
potential 

t/(m)  =  —  6(mm)  :  mm,  (2) 

where  the  tensor  product  mm  and  tensor  double  contraction  A  :  B  are  defined  as 


mimi 

miW2 

miTO3 

mm  = 

rri2mi 

m2m2 

m2m3 

m^rni 

m3m2 

m3m3_ 

A  :  B  =  trace(AB). 

In  (2)  the  constant  b  is  proportional  to  the  normalized  polymer  concentration  and  it  describes 
the  strength  of  inter-molecular  interactions,  and  (mm)  denotes  the  second  moment  of  the 
orientation  distribution: 

(mm)  =  /  mm  p( m)  dm, 

d||m||=l 

where  p( m)  is  the  orientational  probability  density  function  of  the  ensemble  of  polymer  rods, 
i.e. ,  the  probability  density  that  a  polymer  rod  has  direction  m. 
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The  Maier-Saupe  potential  in  a  properly  selected  coordinate  system: 

Let  us  select  a  coordinate  system  such  that  the  second  moment  (mm)  is  diagonal: 

Si  +  g  0  0 

0  S2  +  —  0  , 

o  o  S3  + 1 

where  sj  =  (m?)  —  |  satisfies 


(mm)  = 


(to 


0 


0 

0  (m%)  0 

0  0  (m\) 


Si  +  s2  +  s3  =  (m\  +  to2  +  ml)  -1  =  0. 


In  this  coordinate  system,  the  Maier-Saupe  potential  takes  the  following  expression: 

U(m)  =  —b((m\)m\  +  (m2)rn2  +  (to3)to3) 

= -  ^(sim?  +  s2ml  +  s3m|).  (3) 

The  constant  (—6/3)  does  not  affect  the  dynamics  or  the  stability  of  the  nematic  polymer 
ensemble.  For  simplicity,  we  drop  the  constant  (—6/3).  From  S1+S2+S3  =  0  and  mf+ml+ml  = 
1,  we  have 


S2  -  Si 


2 


S2  ~  Sj 
2 


1  2  1  2 

m1  +  m2  =  1  —  m3. 


Using  these  equations,  we  can  rewrite  the  Maier-Saupe  potential  (3)  as 


where 


U(  m) 


ml 


s3ml 


-6 


(S2  -  Si)  ^2 


(ml  -  m{) 


~Vi 


V2 (ml  -  ml), 


m 


V2 


=  b 


(s2  -  Sj) 
2 


(ra»>  -  5)  • 

h\(ml  -ml). 


(4) 

(5) 

(6) 


Note  that  for  any  orientational  distribution  p(m),  by  selecting  a  proper  coordinate  system,  we 
can  always  express  the  Maier-Saupe  potential  XJ( m)  in  the  form  of  (4)  with  1)  1  and  ?;2  given  by 
(5)  and  (6).  In  particular,  p( m)  does  not  have  to  be  an  equilibrium  distribution. 

The  Boltzmann  distribution: 
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In  the  absence  of  external  field,  the  equilibrium  distributions  of  a  nematic  polymer  ensemble 
are  described  by  the  Boltzmann  distribution  [6] : 


p(BZ)  (m;  m  m)  =  exrt~U^ 

fs  exp [-U (m)]  dm 


exp  [771777,3  +  772(777-2  —  mf) 
fs  exp  [771777.3  +  rj2{m\  —  mf )]  dm 


(7) 


where  S  represents  the  unit  sphere.  It  is  clear  that  an  equilibrium  distribution  is  completely 
specified  by  (771,772).  To  find  an  equilibrium  distribution,  we  need  to  determine  (771,772).  The 
governing  equation  for  (771, 772)  is  obtained  by  combining  the  definition  of  (771, 772)  and  the  Boltz¬ 
mann  distribution  (7): 


77i  =  b\  js  (™3  “  ^  P(BZ) (m; t?i, 772)  dm,, 

V2  =  b\  I  (m2  “  TOi)  P(SZ)(m;  771,772)  dm.  (8) 

J  s 

Note  that  while  the  Maier-Saupe  potential  (4)  is  always  true,  the  Boltzmann  distribution  (7)  is 
valid  only  for  equilibrium  distributions.  It  is  very  important  to  notice  this  difference  between 
(4)  and  (7)  especially  when  we  consider  the  free  energy  of  orientational  distributions  in  the 
stability  discussion. 

Equilibrium  states  of  a  nematic  polymer  ensemble: 

It  has  been  shown  that  in  the  absence  of  external  field,  all  equilibrium  states  of  nematic 
polymers  are  axisymmetric  [8,  24,  35].  Without  loss  of  generality,  we  assume  the  777.3-axis  is  the 
axis  of  symmetry.  The  axisymmetry  implies  Si  =  S2  and  consequently  we  have  772  =  0.  Thus, 
an  equilibrium  state  of  nematic  polymers  is  completely  specified  by  771. 

To  facilitate  the  discussion,  we  use  r(b)  to  denote  the  value  of  771  in  the  equilibrium  state(s) 
corresponding  to  parameter  b.  More  specifically,  we  adopt  the  convention  that  77,  771  and  772 
denote  general  variables  while  r(b)  denotes  specific  value  of  771  in  the  equilibrium  state(s). 

By  definition  r{b )  is  the  solution  of  equation  (8).  With  772  =  0  and  771  simply  denoted  by 
77,  equation  (8)  becomes 


=  fc3  fs  (777.3  ~  |)  exp  (77777!)  dm 
2  fs  exp  (77  777§)  dm 

It  is  now  convenient  to  introduce  the  function  /(??)  defined  by 

=  3  fg  (r»§  -  l)exp(ymj)  dm 
2r/  fs  exp  (77  Tu.g )  c^m 

Then  it  follows  at  once  that  equation  (9)  can  be  written  as  [35] 

\  -  /(» 7))  =  0. 


(9) 


(10) 


(11) 


First,  observe  that  77  =  ?'o(b)  =  0  is  always  a  solution  of  equation  (11),  which  corresponds  to 
U{ m)  =  constant  (independent  of  m),  the  isotropic  equilibrium.  Next,  anisotropic  equilibrium 
states  are  the  solutions  of  the  equation 


f(v)  = 


1 

b’ 


(12) 
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Using  the  spherical  coordinates, 


mi  =  sin  <j>  cos  9 , 
ui2  =  sin</)sin0, 

7713  =  cos  </>,  (13) 


followed  by  the  substitution  u  =  cos  cj>  and  then  applying  integration  by  parts,  we  write  function 

f(v)  as 


3  fg  (cos2  <j)  -  |)  exp  (tj  cos2  cj>)  sin0d<(> 
^  2r]  f*  exp  (77  cos2  (j>)  sin  <j)  d <j> 

3  f- i  (y2  —  g)  exp  (77  u2)  dit 
2^7  f*  exp  (77  it2)  du 
(u2  —  it4)  exp  (77  u2')  dit 
exp  (77  u2)  dit 

fs  (to3  ~  mf)  exP  (v  mi )  dm 
fs  exp  ( 77  m §)  dm 


(14) 


In  [35],  we  have  shown  that  the  function  f(rf)  expressed  in  (14)  satisfies  the  following  properties: 

i-  m  =  h 

2.  lim  /(?7)  =  0  and  lim  f{if)  =  0;  in  other  words,  f(r])  tends  to  zero  as  I77I  goes  to 

00  r] — *■— 00 

infinity. 

3.  /(?7)  attains  its  maximum  at  77*  =  2.1782879748  >  0  where  the  maximum  is  f{r]*)  = 
0.14855559992254  >  0; 

4.  f(ji)  >  0  for  77  <  77*  and  f{if)  <  0  for  77  >  77*. 

The  graph  of  the  function  f(jq)  is  depicted  in  Figure  1. 


Fig.  1  Graph  of  the  function  f(r/).  A  key  feature  of  the  function  /(??)  is  that  f(r/)  is  strictly 
increasing  for  77  <  77*  (i.e.  f{rj)  >  0)  and  strictly  decreasing  for  77  >  77*  (i.e.,  f{rj)  <  0) 


From  these  properties  of  the  function  /(? 7),  we  draw  conclusions  listed  below  related  to 
equation  (12): 

•  b*  =  jjjprj  =  6.7314863965  is  a  critical  value  for  parameter  b. 
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•  For  b  <  b *,  equation  (12)  has  no  solution;  that  is,  there  is  no  anisotropic  state  for  b  <  b*. 

•  At  b  =  6*,  equation  (12)  has  one  solution  r(b*)  =  rf  >  0. 

•  For  b*  <  b  <  4^,  equation  (12)  attains  two  solutions:  rj/(6)  >  rj*  >  0  and  0  <  tm(&)  < 
r/*.  The  equilibrium  states  corresponding  to  positive  order  parameters  rjj(b )  and  t'm  (b)  are 
called  prolate  states. 

•  At  b  =  4y,  equation  (12)  possesses  two  solutions:  >  rj*  >  0  and  =  0. 

•  For  b  >  4y,  equation  (12)  has  two  solutions:  rjj{b)  >  rj*  >  0  and  rz,(6)  <  0.  The 

equilibrium  state  with  negative  order  parameter  r^{b)  is  called  oblate  state. 

Here  we  use  the  subscript  “U”  to  refer  to  the  “Upper”  part  of  the  phase  diagram  where  r  >  rj*, 
the  subscript  “M”  to  represent  the  “Middle”  part  of  the  phase  diagram  where  0  <  r  <  p* ,  and 
the  subscript  “L”  to  denote  the  “Lower”  part  of  the  phase  diagram  where  r  <  0.  The  phase 
diagram  for  nematic  polymers  is  shown  in  Figure  2. 


In  the  above,  we  have  used  b  as  the  independent  variable  and  treated  r  as  a  function  of 
b.  However,  r(b)  is  a  multi-valued  function  of  b  and  for  b  <  b*  the  function  r(b)  is  not  even 
defined.  When  we  study  the  relation  between  b  and  r,  it  is  mathematically  more  convenient 
if  we  use  r  as  the  independent  variable  and  treat  b  as  a  function  of  r  instead.  Then,  b(r)  is  a 
single-valued  function  of  r  and  is  defined  for  all  values  of  r  in  (— oo,  +oo).  The  function  b{r) 
can  be  easily  determined  from  equation  (12)  as 

b (r)  =  -7T-T-  (15) 

f(r) 

Once  we  know  the  function  b(r),  the  branch  rjj{b)  is  simply  the  inverse  function  of  b(r)  for 
r  >  r]*;  the  branch  tm(5)  is  the  inverse  function  of  b(r)  for  0  <r  <rf\  and  the  branch  ri,(b)  is 
the  inverse  function  of  b{r)  for  r  <  0. 

Free  energy  of  an  orientational  distribution: 

The  free  energy  of  the  orientational  distribution  p(m)  is 


G([p})  = 


is  L 


logp(m)  +  -U{ m,  [p]) 


p(m)dm. 


(16) 


Recall  that  for  any  p(m),  by  selecting  a  proper  coordinate  system,  we  can  always  write  the 
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Maier-Saupe  potential  U{ m,  [p])  as 


U{m,  [p])  =  —771 


772(7712  -  ml), 


where  771  =  b |  ((m§)  —  |)  and  772  =  b\{m\  —  ml).  But  in  general,  p( m)  is  not  the  same  as  the 
Boltzmann  distribution 


P{ m)  ^  p^BZ)  (m;  771, 772) 


exp  ??1 777-3  +  r?2(77x|  —  ml) 
fs  exp  [771 777|  +  772(7712  —  mT)]  dm 


It  is  important  to  point  out  that  even  if  we  restrict  our  consideration  to  probability  density 
p(m)  that  has  the  form  of 


,  ,  _  exp  [ <71777.3  +  92(777!  -  ml) _ 

P  fs  exp  [771 777.3  +  92(7712  -  777-1)]  dm ' 

we  still  have  ( 91,92 )  7^  (771,772)  unless  p(m)  is  an  equilibrium  distribution,  which  means  772  =  0 
and  771  =  r(6).  It  is  tempting  to  replace  p(m)  by  p^BZ\ m;  777, 772)  in  (16)  and  formally  define  a 
function  of  (771,772) 


F{  777,772)  =  f 

Js 


logp(BZ)(m;  771,772) 


771  777-3  - 


-  772(77x2  -  777-1) 


Xp(sz)(m.  ^  ^2)dm. 


Unfortunately,  the  assertion  F(rji,  772)  =  G([p])  is  valid  only  when  p(m)  is  an  equilibrium 
distribution.  In  other  words,  at  a  fixed  value  of  b,  the  assertion  #(771,772)  =  G([p])  is  true 
only  at  a  few  (at  most,  three)  equilibrium  states  that  are  well  separated  from  each  other. 
Therefore,  while  we  can  try  to  analyze  the  stability  of  #(771,772)  with  respect  to  perturbations 
in  (771, 772),  the  stability  of  #(?7i,  772)  does  not  inform  us  about  the  stability  of  G([p])  with  respect 
to  perturbations  in  p(m). 


3  Stability  Analysis 

In  this  section  we  present  a  detailed  stability  analysis  on  the  equilibria. 

To  make  our  presentation  easy  to  follow,  we  divide  our  approach  into  five  steps.  In  the  first 
step,  we  minimize  the  free  energy  functional  G([p])  over  a  set  B  of  functions  with  the  second 
moments  fixed,  specified  by  two  parameters  (771,772).  The  resulting  free  energy  is  a  function  of 
two  variables  G (771, 772).  The  stability  of  G([p])  is  equivalent  to  the  stability  of  G (771, 772).  In  the 
second  step,  we  show  that  the  minimum  of  G([p\)  is  attained  at  a  probability  density  p*  of  the 
Boltzmann  form.  The  coefficients  (cji ,  772)  in  the  exponent  of  p*  define  a  one-to-one  mapping 
between  (771,772)  and  (91,92)-  Thus,  we  can  rewrite  G(77i ,  772)  as  ^(91,92)-  The  important  point 
to  note  here  is  that  (91,92)  is  in  general  different  from  (771,772)  unless  p*  is  an  equilibrium 
state,  which  is  an  improper  assumption  in  stability  analysis.  In  the  third  step,  we  compute  the 
Hessian  matrix  of  #(91,92)  with  respect  to  (91,92),  which  is  fortunately  a  diagonal  matrix.  In 
the  fourth  step,  we  consider  the  stability  of  the  isotropic  equilibrium  state.  At  the  last  step  we 
focus  our  attention  on  the  stability  of  the  anisotropic  equilibrium  states. 
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Step  1  We  consider  the  constrained  minimum  over  the  set 

B{vi,V 2)  =  jp(m)  ^(m§)[p]  -  mi)[p]  =  j , 

where  the  probability  density  p( m)  is  involved  in  the  average  in  the  following  way: 

(/(m))[p]  =  /  /(m)p(m)dm. 

Js 

The  minimum  of  G([p])  over  the  set  5(771,772)  is  a  function  of  (771,772): 

£(771,772)=  min  G([p]).  (17) 

p(m)e  #(771 ,7?2) 

From  the  definition  of  G(?7i,  772),  we  see  that  the  stability  of  G ([/>])  with  respect  to  p  is  the  same 
as  the  stability  of  £(771,772)  with  respect  to  (771,772). 

Step  2  We  study  where  G([p])  attains  the  minimum  over  the  set  5(771,772)  and  use  the 
result  to  calculate  £(771,772).  For  p(m)  £  5(771,772),  we  have 


U{ m,  [p])  =  -7?i  |j773  -  -J  ~  772(7712  -  777i), 

■  =  (  (-J?,  (ml  -  -rn{ml-m\)^ 


=  *  +  l4)  ' 


That  is,  the  second  term  in  the  free  energy  G([p])  given  in  (16)  is  actually  a  constant  over  the 
set  5(771,772).  It  follows  that  we  only  need  to  look  at  the  first  term  of  G([p])  in  the  constrained 
minimization  problem.  Let  us  denote  the  first  term  by 


Gi([p])  =  /  log p(m)p(m)dm. 
Js 


The  constrained  minimization  problem  then  becomes 

argmin  Gi  ([/?]). 

p(m)eS(j;i,r/2) 

Suppose  the  constrained  minimum  is  attained  at  p*  (m) .  We  consider  the  perturbed  probability 
density  p(m)  =  p*( m)  +  s  Ap(m)  where  Ap(m)  satisfies  the  following  properties: 

/  Ap(m)dm  =  0, 

Js 

J  (m^  “  A^m)dm  =  °> 

/  (m2  —  mf)Ap(m)dm  =  0.  (19) 

Js 

Gi([p*  +  s  Ap])  attaining  the  minimum  at  s  =  0  implies  dGldp,+a Ap^  =  0,  which  gives  us 

s  s=0 


/  logp*(m)A/3(m)dm  =  0. 

Js 


(20) 
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Since  equation  (20)  is  true  for  all  Ap(m)  satisfying  (19),  it  implies  that 

log(  argmin  Gi([p]))  =  logp*(m)  =  q0  +  91  ( ml  -  ^  )  +  92(m2  -  mf).  (21) 
'  p(m)G-B(?7i,r;2)  '  V  'V 

The  result  of  the  constrained  minimization  given  in  (21)  defines  a  mapping  from  (771,772)  to 
(<7i ,  <72).  The  mapping  from  (771,772)  to  (<71,92)  is  not  expressed  in  an  explicit  form  and  it  is 
not  immediately  obvious  if  the  mapping  is  single- valued.  However,  the  inverse  mapping  from 
(91,92)  back  to  (771,772)  is  single- valued  and  has  a  simple  explicit  expression: 

=  bl  (m>  -  5)  • 

772  =  b^(ml  -  ml),  (22) 

where  the  average  is  with  respect  to  the  probability  density  p*{ m,  <71, 92)  given  by 

p*(m,9i,92)  ee  y  exp  9i  +  q2{ml  -  m\ )  ,  (23) 


Z{qi,  92)  eeee  J  exp  91  ^m2  -  ^  +  92(7712  -  m\)  dm.  (24) 

To  prepare  for  the  computation  of  the  Hessian  matrix,  we  prove  the  following  lemma, 
which  asserts  that  the  mapping  from  (771,772)  to  (91,92)  is  also  single- valued. 

Lemma  1  For  any  fixed  (771,772),  there  is  only  one  set  of  (91,92)  satisfying  (22). 

Proof  We  prove  by  contradiction.  To  accomplish  this,  suppose  (9^ ,  9^ )  and  (9^ ,  ) 

are  different  from  each  other,  and  both  satisfy  (22). 

Let  p*( m)  =  p*(m,  9i1},  92^)  and  Ap(m)  =  p*(m,  q[2] \  q£ ’)  -  p*(m,  q[r \  q£]). 

We  calculate  the  first  and  second  derivatives  of  the  function  G\(\p*  +  s  Ap])  with  respect  to  s 
as  follows: 


dGi([p*  +  s  A p\) 
cl  s 

d2Gi([p*  +  s  Ap]) 


=  /  log(p*  +  s  Ap)Apdm, 

Js 


Is  (p*  +  sAp) 


(Ap)2dm  >  0. 


Since  both  (q^\  q^)  and  (9^,  9^)  satisfy  (22),  Ap(m)  must  satisfy  (19),  which  means 


dGi([p*  +  s  A p\) 
As 


=  f  log(p*)Apdm  =  0. 

Js 


Taylor’s  theorem  tells  us 


Gi{[p*  +  Ap])  =  Gi([p*])  +  1  d2Gl([Pd*g2+  SAP])  _  >  Gi([p*]), 


which  is 


Gi([p*(m,  q[2\  922'*)])  >  ,q^)]). 

Repeating  the  argument  with  the  roles  of  (9^,92^)  and  (q^\q^)  swapped,  we  get 

GiG/tm^^)])  >  Gi([p*(m,9i2),9^2))]), 


(26) 
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which  contradicts  with  (25).  Thus  the  claim  is  established:  {q[1\qi1'1)  and  (q^\q^)  must  be 
the  same.  □ 

We  conclude  that  the  mapping  between  (771,772)  and  (91,92)  is  one-to-one.  As  a  result,  we 
can  express  the  minimum  free  energy  over  the  set  5(771,772)  as  a  function  of  (91,92): 

5(91,92)=  min  G([p]) 
p(m)eB  (771,772) 

=  91  (ml  -0+92  (ml  -  m\)  -  log(Z(9i,  92))  -  ^  Q?7i  +  V2 


0i 771  +  92772^  -  log(Z(9i,  92))  -  i  Q?7i  +  rfe 


(27) 


where  the  mapping  from  (91,92)  to  (771,772)  is  given  in  (22)  and  equations  (16),  (23)  and  (24) 
have  been  invoked. 

The  stability  of  G([p])  with  respect  to  /o(m)  is  the  same  as  the  stability  of  5(91,92)  with 
respect  to  (91,92),  which  is  determined  by  its  Hessian  matrix.  The  Hessian  matrix  (or  simply 
the  Hessian)  is  just  the  square  matrix  of  second-order  partial  derivatives  of  a  function.  It  can 
be  used  to  determine  the  stability  or  instability  of  a  system. 

Step  3  We  now  calculate  the  Hessian  matrix  of  5(91,92).  To  do  so,  we  begin  by 
differentiating  771,  772  and  log(Z(9i,  92))  with  respect  to  91  and  92.  Referring  to  (22),  (23)  and 
(24)  and  by  direct  computation,  we  obtain  the  following  partial  derivatives: 


<9p*(m,  91,  g2) 
dqi 


1 


1 


dp*( m,  91,92) 


dq2 


=  -3  J  75*(m,  9i,  92)  -  (m3  -  -Jp*(m,q1,q2) 

=  (ml  -  (r7i|))p*(m,  9i,  92), 

=  ((7772  “  ml)  -  (ml  -  m\))  p*(xn,q1,q2), 


=  ((m3>  ^  (to3>2)  >  0, 

^  =  b\  i(ml(ml  ~  TOi)>  -  (ml)(m22  -  ml))  , 
=  0  {(ml(ml  -  ml))  -  (ml)  (ml  -  m\))  , 
=  0  (((ml  -  in.?)2)  -  (ml  -  mj)2)  >  0, 


dlog(Z(qllq2)) 

dqi 

d\og(Z(qi,q2)) 


=  '  3  _  3  /  =  36 


/_2  _2\  z 
=  (m2-m1)  =  -772. 


992  v  "  & 

Applying  these  results  to  differentiate  5(91,92)  given  in  (27),  we  get 

05(91,92)  =  2 
“  b 


dqi 

dH(qi,  92)  _  2 
992  b 


1/  ,  9t7i  ,  9t?2 


1/  ,  9t7i  ,  9772 
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To  determine  the  stability  of  an  equilibrium  state,  we  only  need  to  calculate  the  Hessian  matrix 
of  #(51,52)  evaluated  at  the  equilibrium  state. 

At  an  equilibrium  state,  (?7i ,  772) leg  satisfies  equation  (8).  The  mapping  between  (51,52) 
and  (771 , 772)  is  described  in  (22).  Comparing  (8)  and  (22),  we  conclude 

(<7i,92)U  =  (1  n,m)\eq- 

In  addition,  all  equilibrium  states  are  axisynnnetric,  which  amounts  to 

52  | eq  =  0. 

Thus,  at  an  equilibrium  state,  we  have 


(ml  —  ml)  I  =  0, 


(ml(ml  ~m\))\eq  =  0, 


which  leads  to 


Calculating  the  Hessian  matrix  of  H (51, 52)  and  evaluating  it  at  the  equilibrium  state,  we  obtain 

d2H(qi,  52)  2  /  drji\  drn 

dqf  e,“3 b{  dqi  )  dqi  eq  ’ 

d2H(qi,q2)  = 

^51  ^52  eq 

d2H(qi,  52)  I  2  /  <9772  \  dri2  I 


m  I  eq  b\  dq2J  dq2\eq 

Note  that  both  and  are  positive.  Consequently,  the  stability  of  H(qi1q2)  is  determined 
by  the  two  quantities  below: 


Hn  =  1  -  =1~b2  (<TO3>  -  (m3>2) 


Hoo  =  1  - 


=  1  -b-((m2-m2if) 


Since  52|eg  =  0,  we  can  write  ((ml  —  77ia)2)|e<?  as 

/frr.2  „2n2vi  _  Jg(m2  -  in2)2  exp  (gim2)  dm 

'(  2  mi)  fg  exp  (5! m2)  dm 

fj( sin2  (j>)2  f27r( cos2  6  —  sin2  9)2d9^  exp(5i  cos2  (j>)  sin</>d^> 

fo  exp(5i  cos2  4>)  sin  (j)d(j> 

=  ^((1  -ml)2)  . 
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Thus,  we  write  H22  as 


H22  =  l-b-((l-m2)2) 


Step  4  Now  we  turn  our  attention  to  the  stability  of  the  isotropic  state,  which  is 
described  by  (<?i ,  92) |eq  =  (0,0).  With  the  help  of  the  spherical  coordinates  (13),  it  is  straight¬ 
forward  to  derive 


*'3/  \eq 


'-3/1  eq 


Jq  cos2  (j>  sin  cj)  d<j> 
Jq  sin  <J)  d  <j) 


=  L  “2d“  =  L 


74dit  = 


1 

5’ 


((1  -  ml)2)\eq  =  j  (1  -  u2)2du  = 


H“  =  1-b  15’ 

H22  1  b  —  . 

It  follows  readily  that  the  Hessian  matrix  evaluated  at  the  isotropic  equilibrium  state  is  positive 
definite  if  b  <  ^  and  negative  definite  if  b  >  -y.  Therefore,  the  isotropic  equilibrium  state  is 
stable  for  b  <  whereas  it  is  unstable  for  b  >  -y. 

Step  5  Finally,  we  close  our  analysis  by  addressing  the  stability  of  anisotropic  equilib¬ 
rium  states.  An  anisotropic  equilibrium  state  is  described  by  q\  =  r(b)  and  <72  =  0.  First  we 
look  at  Hu.  We  use  the  definition  of  function  f(rj)  given  in  (10)  to  write 


9/(9)  =  2 


3  fs  (m2  -  |)  exp  (1 qm\ )  dm  3 


fs  exp  ( qm §)  dm 


2  1 

=  -  (  rrio  —  - 

2  \  3  3 


(28) 


(gi,qz)=(q,o) 


Note  that  (28)  is  true  for  all  values  of  q,  not  just  at  the  equilibrium  state  q  =  r(b).  Differentiating 
(28)  with  respect  to  q  and  then  evaluating  it  at  the  equilibrium  state,  we  have 


,/(r(6))  +  r(b)f(r(b ))  =  -  ((to£)  -  (ml)2) 


(29) 


eq 


Since  r(b)  is  the  solution  of  eciuation  (12),  it  satisfies  /(?’(&))  =  1/6.  Substituting  this  into  (29) 
yields 

Hn  =  1  -  6  |  ((mi)  -  (m2)2)  =  -br(b)f'(r(b)). 

eq 

Recall  that  in  [35]  we  have  derived  that 

f'(rj)  <  0  for  rj  >  r/*  and  f'(r})  >  0  for  77  <  77* 
and  recall  the  classification  of  three  anisotropic  branches 

ru(b)  >  77*,  0  <  rM(b)  <  V*,  and  rL(b)  <  0, 

we  see  that 


Hn  = 


>  0  for  ru(b), 
<0  for  rM(b), 

>  0  for  r/,(6). 


(30) 
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Next  we  examine  Hii-  Multiplying  H22  by  4/6,  using  1/6  =  f(r(b))  and  the  expression  of  f(rj) 
given  in  (14),  we  have 

\H 22  =  ~  ^  -  m3)2)l(9l,92)  =  (r(6),0) 

=  4(m§  -TO3>|(g1)?2)=(r.(6)i0)  -  ((!  -  m3)2)|(gi,92)=(r(fc)io) 

=  <(1  -m§)(5m§  -  l))|(gii92)=(r(6)i0)  =  h{r(b)), 

where  the  function  6(g)  is  defined  as 


%)  =  <(1  -^3)(5m3  -  !)) 


(9l.92)  =  (q,0) 


(31) 


To  determine  the  sign  of  6(r(6)),  we  study  the  behavior  of  function  6(g).  We  start  by  verifying 

6(0)=  f  (1  —  u2)(5u2  —  l)du  =  0. 

Jo 

Next  we  differentiate  6(g): 

h\q)  =  <(1  -  -  l)(m§  -  (m3»>|(gii,2)=(9,0) 


=  l  ((!  -  ma)(5 ml  -  1)(5 m§  -  1  +  1  -  5(m§)))| (gi>ga)^(g>0) 

=  \  ((!  -  to3)(5to3  -  ®)=(g, o)  +  it1  -  5<to3»%). 


(32) 


Note  that  the  first  term  in  h'(q)  is  always  positive  for  all  values  of  g: 

<(l-mi)(5mi-l)2>|(gii9a)=(gi0)  >0. 

As  a  result,  the  function  6(g)  has  the  property  that  6'(g)  >  0  wherever  6(g)  =  0.  This  property 
along  with  6(0)  =  0  implies 

6(g)  >  0  for  g  >  0  and  6(g)  <  0  for  g  <  0. 


Thus,  for  H22,  we  have 

f  >  0  for  r,j{h), 

H22  =  ^6(r(6))  =  |  >  0  for  rM(b),  (33) 

(  <  0  for  rL(b). 

Putting  the  results  from  (30)  and  (33)  together,  we  arrive  at 

f  (+,+)  for  ru(b), 

(#11,-6622)  =  j  (-,+)  for  rM(b), 

[(+,-)  for  rL (6). 


Therefore,  only  the  highly  aligned  prolate  branch  rjj  (6)  is  stable  while  the  weakly  aligned  prolate 
branch  r^(6)  and  the  oblate  branch  t\l(6)  are  unstable.  More  specifically,  the  branch  ru(b)  is 
stable  to  all  perturbations;  the  branch  Tm(6)  is  unstable  with  respect  to  perturbation  in  gi, 
i.e. ,  axisymmetric  perturbation  along  the  axis  of  symmetry.  Such  perturbation  will  drive  the 
system  either  to  highly  aligned  prolate  branch  rj/(6)  or  the  isotropic  branch.  In  contrast,  the 
oblate  branch  r^(6)  is  unstable  with  respect  to  perturbation  in  q2,  i.e.,  perturbation  away  from 
axisymmetry.  This  kind  of  perturbation  will  carry  the  system  to  the  prolate  state. 
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4  Conclusions  and  Perspective  Views 

The  stability  analysis  performed  in  this  study  confirms  the  well-known  hysteresis  phe¬ 
nomenon  of  the  equilibria  of  the  nematic  liquid  crystalline  polymers  governed  by  the  Smolu- 
chowski  equation  with  the  Maier-Saupe  intermolecular  potential.  In  particular,  the  isotropic 
equilibrium  state  is  stable  for  b  <  15/2  and  it  becomes  unstable  for  b  >  15/2  where  b  is 
proportional  to  the  normalized  polymer  concentration  and  it  describes  the  strength  of  the  in¬ 
termolecular  interactions.  As  for  the  anisotropic  equilibrium  states,  only  the  highly  aligned 
prolate  state  is  stable  while  the  less  aligned  prolate  state  (b*  <  b  <  15/2)  and  the  oblate  state 
( b  >  15/2)  are  unstable.  The  anisotropic  states  occur  when  b  >  b*  where  the  critical  point 
b*  =  6.7314863965. 

While  the  mathematical  study  on  the  Smoluchowski  equation  with  the  Maier-Saupe  po¬ 
tential  is  getting  well-developed,  rigorous  proofs  on  the  similiar  results  for  the  Smoluchowski 
equation  with  the  Onsager  potential  remain  to  be  explored.  We  end  this  paper  by  adopting  a 
perspective  view  provided  by  Donald,  Windle  and  Hanna  [5]:  “It  is  clear  that  ‘self-assembly’, 
where  molecules  are  designed  so  that  they  organise  themselves  into  larger  scale  structures  in 
order  to  achieve  special  properties,  will  be  a  significant  objective  in  materials  science  in  the 
twenty-first  century.  The  demands  of  nanotechnology  require  an  ever  finer  control  over  the 
molecular  arrangements  within  new  materials.  The  fact  that  liquid  crystallinity  itself  is  a  form 
of  orientational  self-assembly,  coupled  with  the  fact  that  the  molecules  in  a  mesophase  can  be 
steered  by  external  fields,  means  that  the  principles  undering  the  science  of  LCPs  can  only  grow 
in  significance  in  the  years  ahead.” 
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